import time
import math  
import random
import numpy as np
import pandas as pd  
import matplotlib.pyplot as plt
  
# 读取Excel文件  
file_path = '机场_单位转换后.xlsx'
df = pd.read_excel(file_path)  
  
# 提取东经和北纬列，转换为列表并合并为一个坐标列表  
coords = df[['东经', '北纬']].values.tolist()  
  
# 将列表转换为NumPy数组  
Point_Map = np.array(coords)
 
DNA_SIZE = len(Point_Map)  # 编码长度（返回行的个数)
POP_SIZE = 1000  # 种群大小
CROSS_RATE = 0.85  # 交叉率
MUTA_RATE = 0.15  # 变异率
Iterations = 1000  # 迭代次数
 
def haversine(lat1, lon1, lat2, lon2):  
    # 将经纬度转换为弧度  
    lat1, lon1, lat2, lon2 = map(math.radians, [lat1, lon1, lat2, lon2])  
  
    # Haversine公式  
    dlon = lon2 - lon1   
    dlat = lat2 - lat1   
    a = math.sin(dlat/2)**2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon/2)**2  
    c = 2 * math.asin(math.sqrt(a))   
    r = 6371  # 地球平均半径，单位为公里  
    return c * r  
  
def distance(DNA):  # 根据DNA的路线计算距离  
    total_distance = 0  
    for i in range(DNA_SIZE - 1):  # 遍历所有城市，除了最后一个  
        from_city = Point_Map[DNA[i]]  
        to_city = Point_Map[DNA[i+1]]  
        # 使用Haversine公式计算距离  
        total_distance += haversine(from_city[1], from_city[0], to_city[1], to_city[0])  # 注意纬度在前，经度在后  
    # 加上从最后一个城市返回第一个城市的距离  
    from_city = Point_Map[DNA[-1]]  
    to_city = Point_Map[DNA[0]]  
    total_distance += haversine(from_city[1], from_city[0], to_city[1], to_city[0])  
    return total_distance

# 用于计算欧几里得距离的情况
# def distance(DNA):  # 根据DNA的路线计算距离
#     dis = 0
#     temp = Point_Map[DNA[0]]
#     for i in DNA[1:]:
#         # sqrt(pow(x-x0,2)+pow(y-y0,2))
#         dis = dis + ((Point_Map[i][0] - temp[0]) ** 2 + (Point_Map[i][1] - temp[1]) ** 2) ** 0.5
#         temp = Point_Map[i]
#     return dis + ((temp[0] - Point_Map[DNA[0]][0]) ** 2 + (temp[1] - Point_Map[DNA[0]][1]) ** 2) ** 0.5
 
 
def getfitness(pop):  # 计算种群适应度，这里适应度用距离的倒数表示
    temp = []
    for i in range(len(pop)):
        temp.append(1 / (distance(pop[i])))
    # 减去最小值是为了防止适应度出现负值
    return temp - np.min(temp)
 
 
def select(pop, fitness):  # 根据适应度选择，以赌轮盘的形式，适应度越大的个体被选中的概率越大
    # print(fitness)
    s = fitness.sum()
    # np.random.choice(a,size,replace,p=None)随机抽取样本a,表示范围，replace=True被抽中后仍有机会被再次抽中，p没抽中的概率
    temp = np.random.choice(np.arange(len(pop)), size=POP_SIZE, replace=True, p=(fitness / s))
    p = []
    for i in temp:
        p.append(pop[i])
 
    return p
 
 
def mutation(DNA, MUTA_RATE):  # 进行变异
    # 两点变异
    if np.random.rand() < MUTA_RATE:  # 以MUTA_RATE的概率进行变异
        mutate_point1 = np.random.randint(0, DNA_SIZE)  # 随机产生一个实数，代表要变异基因的位置
        mutate_point2 = np.random.randint(0, DNA_SIZE)  # 随机产生一个实数，代表要变异基因的位置
        while (mutate_point1 == mutate_point2):  # 保证2个所选位置不相等
            mutate_point2 = np.random.randint(0, DNA_SIZE)  #如果相等将mutate_point2重新进行随机生成位置
        DNA[mutate_point1], DNA[mutate_point2] = DNA[mutate_point2], DNA[mutate_point1]  # 2个所选位置进行互换
    
 
def crossmuta(pop, CROSS_RATE):  # 交叉变异
    new_pop = []
    for i in range(len(pop)):  # 遍历种群中的每一个个体，将该个体作为父代
        n = np.random.rand()
        if n >= CROSS_RATE:  # 大于交叉概率时不发生变异，该子代直接进入下一代
            temp = pop[i].copy()
            new_pop.append(temp)  # 直接进行拷贝
 
        if n < CROSS_RATE:  # 小于交叉概率时发生变异
            list1 = pop[i].copy()
            list2 = pop[np.random.randint(POP_SIZE)].copy()  # 选取种群中另一个个体进行交叉(随机选择)
            status = True
            while status:  # 产生2个不相等的节点，中间部分作为交叉段，采用部分匹配交叉(直到k1<k2的时候才会跳出循环)
                k1 = random.randint(0, len(list1) - 1)
                k2 = random.randint(0, len(list2) - 1)
                if k1 < k2:
                    status = False
 
            k11 = k1  # 保存切片起始的下标
 
            # 先对部分片段进行切片，把切片出来的内容进行交换（完全交换)
            fragment1 = list1[k1: k2]
            fragment2 = list2[k1: k2]
 
            list1[k1: k2] = fragment2
            list2[k1: k2] = fragment1
 
            del list1[k1: k2]  # 删除list1中[k1,k2)的内容
            left1 = list1
 
            # 进行部分匹配的交叉
            offspring1 = []#后代
            #对left1中的每一个位置pos遍历
            for pos in left1:
                #检查它是否存在于frag2中
                if pos in fragment2:
                    #从fragment1中找到对应的基因
                    pos = fragment1[fragment2.index(pos)]
                    #直到基因不再fragment2中为止（遍历fragment2，确保每一个基因都和pos不同)
                    while pos in fragment2:
                        pos = fragment1[fragment2.index(pos)]
                    offspring1.append(pos)
                    continue
                #如何pos不存在fragment2中，那么就直接将其添加到新的后代中
                offspring1.append(pos)
            # 插入新片段
            for i in range(0, len(fragment2)):
                offspring1.insert(k11, fragment2[i])
                k11 += 1
            temp = offspring1.copy()
 
            mutation(temp, MUTA_RATE)  # 进行变异
 
            new_pop.append(temp)  # 把部分匹配交叉后形成的合法个体加入到下一代种群
 
    return new_pop
 
 
def print_info(pop):  # 用于输出结果
    fitness = getfitness(pop)
    maxfitness = np.argmax(fitness)  # 得到种群中最大适应度个体的索引
    # 打印结果
    print("最优的基因型：", pop[maxfitness])
    print("最短距离：", distance(pop[maxfitness]))
    # 按最优结果顺序把地图上的点加入到best_map列表中
    best_map = []
    for i in pop[maxfitness]:
        best_map.append(Point_Map[i])
    best_map.append(Point_Map[pop[maxfitness][0]])
    X = np.array((best_map))[:, 0]
    Y = np.array((best_map))[:, 1]
    # 绘制地图以及路线
    plt.figure()
    plt.rcParams['font.sans-serif'] = ['SimHei']
    plt.scatter(X, Y)
    for dot in range(len(X) - 1):
        plt.annotate(pop[maxfitness][dot], xy=(X[dot], Y[dot]), xytext=(X[dot], Y[dot]))
    plt.annotate('start', xy=(X[0], Y[0]), xytext=(X[0] + 1, Y[0]))
    plt.plot(X, Y)
 
 
if __name__ == "__main__":  # 主循环
    # 生成初代种群pop
    pop = []
    list = list(range(DNA_SIZE))  # 生成[0,DNA_SIZE)的列表
    for i in range(POP_SIZE):  # POP_SIZE是指种群大小，在程序中是一个固定的值(打乱POP_SIZE次之后把结果储存到pop列表中
        random.shuffle(list)  # 随机打乱list，进行初始化操作
        l = list.copy()  # 把list中的数据拷贝到l中
        pop.append(l)  # 将l添加到pop列表中
    best_dis = []
 
    # 获取当前时间(算法开始时间)
    start_time = time.time()
 
    # 进行选择，交叉，变异，并把每代的最优个体保存在best_dis中
    for i in range(Iterations):  # 迭代N代
        pop = crossmuta(pop, CROSS_RATE)  # CROSS_RATE交叉率
        fitness = getfitness(pop)  # 得到适应度种群的适应度
        # 更新最差适应度
        # print(np.min(fitness))
        tmpfitness = np.max(fitness)
 
        maxfitness = np.argmax(fitness)  # 返回数值最大的索引
        best_dis.append(distance(pop[maxfitness]))
        pop = select(pop, fitness)  # 选择生成新的种群（适应度最大的)
 
    #averageFitness = sumFitness / sumCount
    # 获取当前时间(算法结束时间)
    end_time = time.time()
 
    print_info(pop)  # 打印信息

